covid-immigrant-foodworkers

Import and format data

Import IPUMS

Use of data from IPUMS USA is subject to conditions including that users should
cite the data appropriately. Use command `ipums_conditions()` for more details.

|=                                                                                                         |   1%    6 MB
|===                                                                                                       |   3%   13 MB
|=====                                                                                                     |   5%   19 MB
|=======                                                                                                   |   7%   26 MB
|=========                                                                                                 |   9%   32 MB
|===========                                                                                               |  10%   39 MB
|=============                                                                                             |  12%   45 MB
|===============                                                                                           |  14%   52 MB
|=================                                                                                         |  16%   58 MB
|===================                                                                                       |  18%   65 MB
|=====================                                                                                     |  19%   71 MB
|=======================                                                                                   |  21%   78 MB
|=========================                                                                                 |  23%   84 MB
|===========================                                                                               |  25%   91 MB
|============================                                                                              |  27%   97 MB
|==============================                                                                            |  28%  104 MB
|================================                                                                          |  30%  110 MB
|==================================                                                                        |  32%  117 MB
|====================================                                                                      |  34%  123 MB
|======================================                                                                    |  36%  130 MB
|========================================                                                                  |  37%  136 MB
|==========================================                                                                |  39%  143 MB
|============================================                                                              |  41%  149 MB
|==============================================                                                            |  43%  156 MB
|================================================                                                          |  45%  162 MB
|==================================================                                                        |  46%  169 MB
|====================================================                                                      |  48%  175 MB
|======================================================                                                    |  50%  182 MB
|=======================================================                                                   |  52%  188 MB
|=========================================================                                                 |  54%  195 MB
|===========================================================                                               |  55%  201 MB
|=============================================================                                             |  57%  208 MB
|===============================================================                                           |  59%  214 MB
|=================================================================                                         |  61%  221 MB
|===================================================================                                       |  63%  227 MB
|=====================================================================                                     |  64%  234 MB
|=======================================================================                                   |  66%  240 MB
|=========================================================================                                 |  68%  247 MB
|===========================================================================                               |  70%  253 MB
|=============================================================================                             |  72%  260 MB
|===============================================================================                           |  73%  266 MB
|=================================================================================                         |  75%  273 MB
|==================================================================================                        |  77%  279 MB
|====================================================================================                      |  79%  286 MB
|======================================================================================                    |  81%  292 MB
|========================================================================================                  |  82%  299 MB
|==========================================================================================                |  84%  305 MB
|============================================================================================              |  86%  312 MB
|==============================================================================================            |  88%  318 MB
|================================================================================================          |  90%  325 MB
|==================================================================================================        |  91%  331 MB
|====================================================================================================      |  93%  338 MB
|======================================================================================================    |  95%  344 MB
|========================================================================================================  |  97%  351 MB
|==========================================================================================================|  99%  357 MB
|===========================================================================================================| 100%  360 MB

PUMA shapefiles

OGR data source with driver: ESRI Shapefile 
Source: "/Users/chriszs/Desktop/work/covid-immigrant-foodworkers/data/ipums_puma_2010/ipums_puma_2010.shp", layer: "ipums_puma_2010"
with 2378 features
It has 7 fields
Integer64 fields read as strings:  GISMATCH 

County shapefiles

OGR data source with driver: ESRI Shapefile 
Source: "/Users/chriszs/Desktop/work/covid-immigrant-foodworkers/data/tl_2019_us_county/tl_2019_us_county.shp", layer: "tl_2019_us_county"
with 3233 features
It has 17 fields
Integer64 fields read as strings:  ALAND AWATER 

COVID county cases and deaths

covid <- read_csv("data/us-counties.csv",
                  col_types = cols(.default = "?", date = "D")) %>%
  rename(county_code = fips) %>%
  # Filter to the latest data for each county
  group_by(county_code) %>%
  slice(which.max(date))

covid

County populations

Read-in Census API key

census_key <- (readLines(key_file)[1])

Create list of all US states to iterate through in the API calls

us <- unique(fips_codes$state)[1:51]

us

 [1] "AL" "AK" "AZ" "AR" "CA" "CO" "CT" "DE" "DC" "FL" "GA" "HI" "ID" "IL" "IN" "IA" "KS" "KY" "LA" "ME" "MD" "MA" "MI"
[24] "MN" "MS" "MO" "MT" "NE" "NV" "NH" "NJ" "NM" "NY" "NC" "ND" "OH" "OK" "OR" "PA" "RI" "SC" "SD" "TN" "TX" "UT" "VT"
[47] "VA" "WA" "WV" "WI" "WY"

Download population data from the 2014-2018 ACS

PUMA-to-county crosswalk

county_to_puma_crosswalk <- read_csv("data/geocorr2018.csv", skip = 1) %>%
  mutate(puma_code = paste0(`State code`, `PUMA (2012)`)) %>%
  select(
    puma_code,
    county_code = `County code (2014)`,
    state = `State abbreviation`,
    county = `2014 county name`,
    puma = `PUMA12 name`,
    allocation_factor = `puma12 to county14 allocation factor`
  )

county_to_puma_crosswalk

Join the COVID, population and crosswalk data frames

crosswalk <- list(covid, population, county_to_puma_crosswalk) %>%
  reduce(full_join, by = "county_code") %>%
  mutate(
    puma = paste(puma, state.x, sep = ", "),
    cases_per_100k = round(cases / population * 100000, digits = 0),
    deaths_per_100k = round(deaths / population * 100000, digits = 0)
  ) %>%
  select(
    1,
    state = state.x,
    county_code,
    county = county.y,
    puma_code,
    puma,
    population,
    cases,
    deaths,
    cases_per_100k,
    deaths_per_100k,
    allocation_factor
  )

crosswalk

Filter to the food production industries we’re interested in

food_production_industries <-
  list("0170", # Crop production
       "0180", # Animal production and aquaculture
       "0290", # Support activities for agriculture and forestry
       "1070", # Animal food, grain and oilseed milling
       "1080", # Sugar and confectionery products
       "1090", # Fruit and vegetable preserving and specialty food manufacturing
       "1170", # Dairy product manufacturing
       "1180", # Animal slaughtering and processing
       "1270", # Bakeries and tortilla manufacturing, except retail bakeries
       "1280", # Seafood and other miscellaneous foods; n.e.c.
       "1290") # Not specified food industries

ipums_food_production <- ipums %>%
  filter(industry_code %in% food_production_industries)

Calculate employment by industry-occupation grouping

ipums_food_production_employment_by_ind_occ <-
  ipums_food_production %>%
  group_by(industry, occupation) %>%
  summarize(num_employees = sum(PERWT)) %>%
  pivot_wider(
    id_cols = c(industry, occupation),
    values_from = num_employees,
    names_prefix = "num_employees"
  ) %>%
  replace(is.na(.), 0) %>%
  # Exclude from our analysis ind-occ groupings
  # with fewer than 2,500 employees nationally
  filter(num_employees >= 2500)

Export the data

write_csv(
  ipums_food_production_employment_by_ind_occ,
  "data/exported/ipums_food_production_employment_by_ind_occ.csv"
)

Filter the industry-occupation pairs

Import data frame of food production industry-occupation pairs

ipums_food_production_employment_by_ind_occ <-
  read_csv("data/ipums_food_production_total_employment_by_ind_occ.csv")

Join the industry-occupation pairs to the ipums data and filter to just the industry-occupation pairs we want

ipums_food_production_selected_occupations <-
  ipums_food_production %>%
  inner_join(ipums_food_production_employment_by_ind_occ,
             by = c("industry", "occupation")) %>%
  filter(exclude == FALSE) %>%
  select(-num_employees, -exclude)

ipums_food_production_selected_occupations

And join with the pumas data frame to get the PUMA names

ipums_food_production_selected_occupations <-
  ipums_food_production_selected_occupations %>%
  left_join(pumas, by = "puma_code") %>%
  mutate(puma = paste(puma, state, sep = ", ")) %>%
  select(1:3,
         13,
         4:12)

ipums_food_production_selected_occupations

citizenship_lookup = tibble(
  citizenship_code = factor(c(
    0, 1, 2, 3, 4, 5
  )),
  citizenship = c(
    "native_citizen", # n/a
    "native_citizen", # born_abroad_citizen
    "naturalized_citizen", # naturalized_citizen
    "non_citizen", # non_citizen
    "non_citizen_with_papers", # non_citizen_with_papers
    "foreign_born_status_unknown"
  )
)

ethnicity_lookup = tibble(
  ethnicity_code = factor(c(
    0, 1, 2, 3, 4, 9
  )),
  ethnicity = c(
    "non_hispanic",
    "hispanic", # mexican
    "hispanic", # puerto_rican
    "hispanic", # cuban
    "hispanic", # other
    "not_reported"
  )
)

ipums_food_production_selected_occupations_mapped <-
  ipums_food_production_selected_occupations %>%
  select(-ethnicity,-citizenship) %>%
  left_join(ethnicity_lookup, by="ethnicity_code") %>%
  left_join(citizenship_lookup, by="citizenship_code")

pct <- function (part,whole) {
  round(part / whole * 100, digits = 0)
}

calc_percents <- function (rows) {
  rows %>%
  mutate(
    total_native_citizen = (
      total_native_citizen_non_hispanic +
      total_native_citizen_hispanic
    ),
    total_naturalized_citizen = (
      total_naturalized_citizen_non_hispanic +
      total_naturalized_citizen_hispanic
    ),
    total_non_citizen = (
      total_non_citizen_non_hispanic +
      total_non_citizen_hispanic
    ),
    total_foreign_non_hispanic = (
      total_naturalized_citizen_non_hispanic +
      total_non_citizen_non_hispanic
    ),
    total_foreign_hispanic = (
      total_naturalized_citizen_hispanic +
      total_non_citizen_hispanic
    ),
    total_foreign = (
      total_foreign_non_hispanic +
      total_foreign_hispanic
    ),
    total_employees = (
        total_native_citizen +
        total_naturalized_citizen +
        total_non_citizen
    ),
    pct_native_citizen_non_hispanic = 
      pct(total_native_citizen_non_hispanic, total_employees),
    pct_native_citizen_hispanic =
      pct(total_native_citizen_hispanic, total_employees),
    pct_total_foreign_non_hispanic =
      pct(total_foreign_non_hispanic, total_employees),
    pct_total_foreign_hispanic =
      pct(total_foreign_hispanic, total_employees),
    pct_naturalized_citizen_non_hispanic =
      pct(total_naturalized_citizen_non_hispanic, total_employees),
    pct_naturalized_citizen_hispanic =
      pct(total_naturalized_citizen_hispanic, total_employees),
    pct_non_citizen_non_hispanic =
      pct(total_non_citizen_non_hispanic, total_employees),
    pct_non_citizen_hispanic =
      pct(total_non_citizen_hispanic, total_employees),
    pct_native_citizen =
      pct(total_native_citizen, total_employees),
    pct_total_foreign =
      pct(total_foreign, total_employees),
    pct_naturalized_citizen =
      pct(total_naturalized_citizen, total_employees),
    pct_non_citizen =
      pct(total_non_citizen, total_employees)
  )
}

ipums_food_production_selected_occupations_mapped

Calculate the national-level citizenship and ethnicity data for these industry-occupation pairs

ipums_food_production_by_industry_occupation <-
  ipums_food_production_selected_occupations_mapped %>%
  group_by(industry_code,
           industry,
           occupation_code,
           occupation,
           citizenship,
           ethnicity) %>%
  summarize(num_employees = sum(PERWT)) %>%
  pivot_wider(
    id_cols = c(industry_code, industry, occupation_code, occupation),
    names_from = c(citizenship, ethnicity),
    values_from = num_employees,
    names_prefix = "total_"
  ) %>%
  replace(is.na(.), 0) %>%
  calc_percents()

ipums_food_production_by_industry_occupation

Export the data

write_csv(
  ipums_food_production_by_industry_occupation,
  "data/exported/ipums_food_production_by_industry_occupation.csv"
)

Group the data by PUMA

ipums_food_production_by_puma <-
  ipums_food_production_selected_occupations_mapped %>%
  group_by(state_code,
           state,
           puma_code,
           puma,
           citizenship,
           ethnicity) %>%
  summarize(num_employees = sum(PERWT, na.rm = TRUE)) %>%
  pivot_wider(
    id_cols = c(state_code, state, puma_code, puma),
    names_from = c(citizenship, ethnicity),
    values_from = num_employees,
    names_prefix = "total_"
  ) %>%
  replace(is.na(.), 0) %>%
  calc_percents()

ipums_food_production_by_puma %>% select(puma_code,total_foreign)

Export the data

write_csv(
  ipums_food_production_by_puma,
  "data/exported/ipums_food_production_by_puma.csv"
)

Group the data by PUMA and industry and occupation

ipums_food_production_by_puma_ind_occ <-
  ipums_food_production_selected_occupations_mapped %>%
  group_by(
    state_code,
    state,
    puma_code,
    puma,
    industry_code,
    industry,
    occupation_code,
    occupation,
    citizenship,
    ethnicity
  ) %>%
  summarize(num_employees = sum(PERWT)) %>%
  pivot_wider(
    id_cols = c(
      state_code,
      state,
      puma_code,
      puma,
      industry_code,
      industry,
      occupation_code,
      occupation
    ),
    names_from = c(citizenship, ethnicity),
    values_from = num_employees,
    names_prefix = "total_"
  ) %>%
  replace(is.na(.), 0) %>%
  calc_percents()

Export the data

write_csv(
  ipums_food_production_by_puma_ind_occ,
  "data/exported/ipums_food_production_by_puma_ind_occ.csv"
)

ipums_food_production_by_puma

Join the employment data to the crosswalk data frame

ipums_food_production_by_puma

ipums_food_production_by_puma_crosswalked_to_county <-
  ipums_food_production_by_puma %>%
  left_join(crosswalk, by = "puma_code") %>%
  mutate(
    total_employees = round(
      total_employees * allocation_factor,
      digits = 0
    ),
    total_native_citizen_non_hispanic = round(total_native_citizen_non_hispanic *
                                          allocation_factor, digits = 0),
    total_native_citizen_hispanic = round(total_native_citizen_hispanic *
                                      allocation_factor, digits = 0),
    total_foreign_non_hispanic = round(total_foreign_non_hispanic *
                                         allocation_factor, digits = 0),
    total_foreign_hispanic = round(total_foreign_hispanic *
                                     allocation_factor, digits = 0),
    total_naturalized_citizen_non_hispanic = round(total_naturalized_citizen_non_hispanic
                                             * allocation_factor, digits = 0),
    total_naturalized_citizen_hispanic = round(total_naturalized_citizen_hispanic
                                         * allocation_factor, digits = 0),
    total_non_citizen_non_hispanic = round(total_non_citizen_non_hispanic
                                     * allocation_factor, digits = 0),
    total_non_citizen_hispanic = round(total_non_citizen_hispanic
                                 * allocation_factor, digits = 0),
    total_native_citizen = total_native_citizen_non_hispanic + total_native_citizen_hispanic,
    total_foreign = total_foreign_non_hispanic + total_foreign_hispanic,
    total_naturalized_citizen = total_naturalized_citizen_non_hispanic +
      total_naturalized_citizen_hispanic,
    total_non_citizen = total_non_citizen_non_hispanic + total_non_citizen_hispanic
  ) %>%
  select(1,
         state = state.x,
         3,
         puma = puma.x,
         32:33,
         35:40,
         5:17)

ipums_food_production_by_puma_crosswalked_to_county

ipums_food_production_by_puma_crosswalked_to_county %>% filter(county_code == "06053") %>% select(allocation_factor, total_employees, total_foreign, total_non_citizen)

Join the employment data to the crosswalk data frame

ipums_food_production_by_puma_ind_occ_crosswalked_to_county <-
  ipums_food_production_by_puma_ind_occ %>%
  left_join(crosswalk, by = "puma_code") %>%
  mutate(
    total_employees = round(
      total_employees * allocation_factor,
      digits = 0
    ),
    total_native_citizen_non_hispanic = round(total_native_citizen_non_hispanic
                                        * allocation_factor, digits = 0),
    total_native_citizen_hispanic = round(total_native_citizen_hispanic
                                    * allocation_factor, digits = 0),
    total_foreign_non_hispanic = round(total_foreign_non_hispanic
                                       * allocation_factor, digits = 0),
    total_foreign_hispanic = round(total_foreign_hispanic
                                   * allocation_factor, digits = 0),
    total_naturalized_citizen_non_hispanic = round(total_naturalized_citizen_non_hispanic
                                             * allocation_factor, digits = 0),
    total_naturalized_citizen_hispanic = round(total_naturalized_citizen_hispanic
                                         * allocation_factor, digits = 0),
    total_non_citizen_non_hispanic = round(total_non_citizen_non_hispanic
                                     * allocation_factor, digits = 0),
    total_non_citizen_hispanic = round(total_non_citizen_hispanic
                                 * allocation_factor, digits = 0),
    total_native_citizen = total_native_citizen_non_hispanic + total_native_citizen_hispanic,
    total_foreign = total_foreign_non_hispanic + total_foreign_hispanic,
    total_naturalized_citizen = total_naturalized_citizen_non_hispanic +
      total_naturalized_citizen_hispanic,
    total_non_citizen = total_non_citizen_non_hispanic + total_non_citizen_hispanic
  ) %>%
  select(1,
         state = state.x,
         3,
         puma = puma.x,
         36:37,
         5:8,
         39:44,
         9:21)

ipums_food_production_by_puma_crosswalked_to_county

Group the data by county

# This step is required to deal with situations
# where multiple PUMAs split a single county
ipums_food_production_by_county <-
  ipums_food_production_by_puma_crosswalked_to_county %>%
  group_by(state_code, state, county_code, county) %>%
  calc_percents()

Export the data

write_csv(
  ipums_food_production_by_county,
  "data/exported/ipums_food_production_by_county.csv"
)

Group the data by industry and occupation and county

ipums_food_production_by_puma_ind_occ_crosswalked_to_county

# This step is required to deal with situations
# where multiple PUMAs split a single county
ipums_food_production_by_county_ind_occ <-
  ipums_food_production_by_puma_ind_occ_crosswalked_to_county %>%
  group_by(
    state_code,
    state,
    county_code,
    county,
    industry_code,
    industry,
    occupation_code,
    occupation
  ) %>%
  calc_percents()

ipums_food_production_by_county_ind_occ

Export the data

write_csv(
  ipums_food_production_by_county_ind_occ,
  "data/exported/ipums_food_production_by_county_ind_occ.csv"
)

Analyze the data

What is the total number of workers employed in these industries?

sum(ipums_food_production_by_puma$total_employees)

What is the total number of foreigners employed in these industries?

sum(ipums_food_production_by_puma$total_foreign)

What’s that as a percent?

sum(ipums_food_production_by_puma$total_foreign) /
  sum(ipums_food_production_by_puma$total_employees) * 100

How many of the jobs have a higher proportion of foreigners employed than the national average of 17.1%?

ipums_food_production_by_industry_occupation %>%
  filter(pct_total_foreign > 17.1) %>%
  select(
    industry,
    occupation,
    total_employees,
    pct_total_foreign,
    pct_total_foreign_hispanic
  ) %>%
  arrange(desc(total_employees))

What is the total number of Hispanic foreigners employed in these industries?

sum(ipums_food_production_by_puma$total_foreign_hispanic)

What’s that as a percent of all foreign workers in these industries?

sum(ipums_food_production_by_puma$total_foreign_hispanic) /
  sum(ipums_food_production_by_puma$total_foreign) * 100

Which PUMAs have the highest number of immigrant foodworkers?

ipums_food_production_by_puma %>%
  arrange(desc(total_foreign)) %>%
  select(
    state,
    puma,
    total_employment_selected_industries,
    total_foreign,
    total_foreign_hispanic
  ) %>%
  head(23) # 2,341 PUMAs / 100 = 23.41

And how many immigrant workers are in these PUMAs?

ipums_food_production_by_puma %>%
  arrange(desc(total_foreign)) %>%
  select(puma,
         total_employment_selected_industries,
         total_foreign,
         total_foreign_hispanic) %>%
  head(23) %>%
  ungroup() %>%
  summarize(sum(total_foreign))

And how many of these foreign-born workers are Hispanic?

ipums_food_production_by_puma %>%
  arrange(desc(total_foreign)) %>%
  select(puma,
         total_employment_selected_industries,
         total_foreign,
         total_foreign_hispanic) %>%
  head(23) %>%
  ungroup() %>%
  summarize(sum(total_foreign_hispanic))

Which counties have the highest number of immigrant foodworkers?

ipums_food_production_by_county %>%
  arrange(desc(total_foreign)) %>%
  select(
    state,
    county,
    total_employment_selected_industries,
    total_foreign,
    total_foreign_hispanic
  ) %>%
  head(31) # 3,142 counties / 100 = 31.42

And how many immigrant workers are in these counties?

ipums_food_production_by_county %>%
  arrange(desc(total_foreign)) %>%
  select(
    county,
    total_employment_selected_industries,
    total_foreign,
    total_foreign_hispanic
  ) %>%
  head(31) %>%
  ungroup() %>%
  summarize(sum(total_foreign))

And how many of these foreign-born workers are Hispanic?

ipums_food_production_by_county %>%
  arrange(desc(total_foreign)) %>%
  select(
    county,
    total_employment_selected_industries,
    total_foreign,
    total_foreign_hispanic
  ) %>%
  head(31) %>%
  ungroup() %>%
  summarize(sum(total_foreign_hispanic))

What is the distribution of foreign workers in these industries by state?

foodworkers_by_state <- ipums_food_production_by_puma %>%
  group_by(state) %>%
  summarize(
    state_total_employment_selected_industries =
      sum(total_employees),
    state_total_foreign = sum(total_foreign),
    state_pct_total_foreign = round(
      sum(total_foreign) / sum(total_employees) * 100
    ),
    2,
    state_total_foreign_hispanic = sum(total_foreign_hispanic),
    state_pct_total_foreign_hispanic = round(
      sum(total_foreign_hispanic) /
        sum(total_employees) * 100
    ),
    2
  ) %>%
  mutate(
    rank_total_employment = dense_rank(desc(
      state_total_employment_selected_industries
    )),
    rank_total_foreign = dense_rank(desc(state_total_foreign)),
    rank_total_foreign_hispanic = dense_rank(desc(state_total_foreign_hispanic))
  ) %>%
  select(
    state,
    state_total_employment_selected_industries,
    state_total_foreign,
    state_pct_total_foreign,
    state_total_foreign_hispanic,
    state_pct_total_foreign_hispanic,
    rank_total_employment,
    rank_total_foreign,
    rank_total_foreign_hispanic
  ) %>%
  arrange(desc(state_total_employment_selected_industries))

foodworkers_by_state

What proportion of total foodworker employment do the top 10 states account for?

foodworkers_by_state %>%
  arrange(desc(state_total_employment_selected_industries)) %>%
  head(10) %>%
  summarize(sum(state_total_employment_selected_industries))

What is the distribution of foreign workers in these industries by state?

foodworkers_by_state %>%
  summarize(sum(state_total_employment_selected_industries))

Where do the “Red Zone” states stand?

red_zone <-
  data.frame(
    "state" = c(
      "Alabama/AL",
      "Arizona/AZ",
      "Arkansas/AR",
      "California/CA",
      "Florida/FL",
      "Georgia/GA",
      "Idaho/ID",
      "Iowa/IA",
      "Kansas/KS",
      "Louisiana/LA",
      "Mississippi/MS",
      "Missouri/MO",
      "Nevada/NV",
      "North Carolina/NC",
      "North Dakota/ND",
      "Oklahoma/OK",
      "South Carolina/SC",
      "Tennessee/TN",
      "Texas/TX",
      "Utah/UT",
      "Wisconsin/WI"
    ),
    "red_zone" = TRUE
  )

red_zone

Join with the red zone states data frame

foodworkers_by_state <- foodworkers_by_state %>%
  left_join(red_zone, by = "state")

Where do the Red Zone states stand?

foodworkers_by_state %>%
  filter(red_zone == TRUE)

Export the states

write_csv(foodworkers_by_state,
          "data/exported/foodworkers_by_state.csv")